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Abstract 

The objectives of this paper is l)to describe an efficient numerical mesh generation 
technique capable of creating single grids around complex winged configurations, 
2)use the above grid generation scheme to simulate the flow field around a realis- 
tic three dimensional configuration by solving the Reynolds averaged Navier-Stokes 
equations using a time dependent implicit procedure. Numerical results are pre- 
sented for the flow around the complete configuration of the shuttle orbiter and a 
canard projectile at supersonic Mach numbers and zero angle of attack. 

Introduction 

The present paper describes a numerical mesh generation technique (Ref. jlj) and 
an implicit finite-difference procedure (Refs. [2,3]) for simulating viscous flow about 
low- aspect-ratio wing body configurations using a single grid strategy. Meshes are 
generated for complex three-dimensional winged configurations using equations of 
mesh line orthogonality and cell volume which form a hyperbolic partial differential 
equation set. The grid generation equations are readily solved in the computational 
domain by marching outward from the body surface. In supersonic flow, the three- 
dimensional Reynolds-averaged Navier Stokes equations are efficiently solved by 
segmenting the computational domain into multiple regions([3]).This is possible 
because the outflow boundary of each segment can be located in supersonic areas 
to avoid upstream influence between two contiguous regions. 

The numerical procedures have been used to compute the viscous high Reynolds 
number flow around a complete configuration of the shuttle orbiter (including the 
canopy, wing. OMS pods and vertical tail) and a canard projectile at supersonic 
Mach numbers and zero angle of attack. In the following sections, we briefly describe 
the grid generation procedure, the flow simulation procedure, and the computed 
results. 

Numerical Grid Generation Procedure 

Although, there are various grid generation techniques available, (see e.g. [4]), 
the hyperbolic grid generation method (jl],[5|) was selected in the present applica- 
tions because it has several desirable features. With the hyperbolic grid generation 



method, the user is able to specify an arbitrary surface point distribution, and the 
numerically generated grid is orthogonal or is as close to being orthogonal as possi- 
ble. The hyperbolic grid generation procedure is also efficient in terms of computer 
time and memory. 

Grid Generation Equations 

Body-fitted curvilinear coordinates of £ = £ (x, y, z), 77 — 77 (x, y, z), and £ = f 
(x, y, z), where f is the outward direction, are sought to simplify the application 
of finite-difference viscous flow field simulation procedures. Following Ref. [l], two 
orthogonality relations and a cell volume constraint are used as grid generation 
equations to provide x, y and z in terms of £, y and f . The resulting equations are: 


ztXf + VtVi + H Z S = 0 (la) 

X**; + Vr,y- + Z V Z( = 0 (lb) 
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The first two equations represent orthogonality relations between £ and f lines and 
between y and f lines, whereas the last equation is the cell volume or the finite 
Jacobian constraint. 

Numerical Solution Algorithm 

The nonlinear system of grid generation equations given by Eqs. (1) are linearized 
about a known state ro, and are solved using a noniterative implicit finite-difference 
scheme. An unconditionally stable implicit scheme is chosen so that the marching 
step size in $ can be arbitrarily selected based only on considerations of accurately 
generating the grid. An iterative solution of the nonlinear grid generation equations 
is avoided by expanding the equations about the previous marching step. As a 
consequence, the resulting equations are solved with the nearby known state 0 
taken from the previous f step. 

Let A£ = Ay = Af = 1 such that £ = j - T, y = k — 1, and f = l — 1. Central 
spatial differencing of Equations (l) in £ and y with two-point backward implicit 
differencing in f leads to 


~ n) + Bi 6 v {fi + 1 - fi) + C;V ? fi + 1 = gi + i (3) 

where / is the identity matrix, and Ai,Bi and Cj are 3x3 coefficient matrices that 
are formed from the local linearization of Eqn. (2). The vector gi + j contains the 
user specified cell volumes 
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( 0 

9l + 1 = 0 

\AV l+1 

The central difference operators used in £ and t] are given by 


h'i - f,> * 2 r, ~ 1 ■ Kn = ft+1 ~ rt -‘ 


whereas the backward differencing used in £ is represented as 


V f rj + i = r t+1 - r t 

To reduce the inversion cost, the difference equations are multiplied through by 
Cf 1 and are approximately factored as 

(/ + Cf'BiSJil + Cf i A l S i )(r l+ j - r t ) = C^g l+ , (4) 

so that f/ + i is obtained by solving sequences of one-dimensional-like block tridiag- 
onal systems 


{J + C/ l BiS T1 )g t + j — C l l gi+i 

(5a) 

{1 + C l 1 Ai6z)¥ <rt + 1 — g /+i 

(5b) 

ri+ 1 = Ft -r V«-r ; + 1 

(5c) 


In practice numerical dissipation terms are added in the £ and r? directions and 
for body shapes that have concavity, super implicitness is built into the algorithm 
following Ref. |6|. 

Initial and Boundary Conditions 

The body surface is mapped to the bottom of the computational cube and the 
grid is obtained by marching outward from the body surface. The surface grid 
can be arbitrarly specifed by the user. For relatively simple body shapes, the 
user can generate an orthogonal grid along the body surface. This process assures 
that the resulting grid will be nearly orthogonal in all directions. For complicated 
configurations, however, the surface grid is generated along cross sections of constant 
axial distance. Depending on the body shape, the £ and r) lines can be far from 
being orthogonal. However, this does not introduce any inconsistency because the 
orthogonality between the £ and r] is not enforced in the grid generation procedure. 

Boundary conditions in the £ and rj directions at each marching step are needed 
to create the grid. Periodicity is always used in the rj (circumferential) direction. A 
stable boundary condition procedure is needed at £ = 0 and £ =£ mai for creating 
“O-O” meshes around closed bodies(Fig. (la)) and “C-0” or “H-0” meshes around 
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open bodies(Fig. (lb, lc)) As shown in Ref. [lj, the eigenvalues of either A or 
Cf l B are of the form (<7,0,-a) where o is real. Consequently, two stable boundary 
procedures have been used: 1) specify one combination of the dependent variables 
and use one sided differencing for the two remaining governing equations at the 
boundary, or 2)specify two combinations of the dependent variables and use one 
sided differencing for a third governing equation at the boundary. The first approach 
is used at a boundary surface (e.g. £ = £ max in Figs, (lb, lc)) whereas the second 
approach is used at an axis (e.g. £ = 0 in Figs, (la, lb)). 

Since, the grid is obtained by marching outward from the body, an exact location 
for the outer boundary can not be prescribed. Nevertheless, an approximate location 
for the outer boundary can be specified through proper selection of the cell volumes 
used in Eqn. (lc). 

Flow Simulation Procedure 

The main objective of the present work is to demonstrate an efficient numerical 
procedure for computing three-dimensional viscous flow over winged configurations 
(e.g., the shuttle orbiter) at supersonic free-stream Mach numbers. Although the 
flow is predominantly supersonic, some sections of the geometry generate local re- 
gions of embedded subsonic flow and/or axial separation. In such regions, the flow 
equations are elliptic, precluding the use of the more efficient space-marching tech- 
niques such as the parabolized Navier-Stokes (PNS) codes (e.g. Ref. 7 ). Instead, 
a time-dependent continuation procedure is used to find the steady-state solution 
of the three-dimensional conservation-law form of the Reynolds-averaged Navier- 
Stokes equations with the thin layer approximation. 

The numerical procedure employs the vectorized Beam- Warming, approximate 
factorization, finite-difference scheme as presented in Refs. [2,3j. The algorithm is 
implicit, non-iterative and uses second or fourth-order accurate central-differencing 
in the spatial directions. A generalized coordinate transformation is employed to 
map the three- dimensional region between the body and the outer boundary into a 
cubical comput ational domain. The grid points are clustered in the direction normal 
to the body surface in such a manner as to adequately resolve both the viscous and 
the inviscid regions. The algebraic model of Ref. [8] is used for turbulent computa- 
tions. The bow shock can either be fitted or captured while internal discontinuities 
are always captured. Calculations are performed using bilateral symmetry with the 
computational domain extending from the windward to the leeward side. 

In order to facilitate the solution of the Navier-Stokes equations around complex 
bodies, a segmentation procedure in which the flow field is divided into multiple 
regions is used. The boundary between two adjacent regions is carefully located in 
areas where the outer inviscid flow is supersonic and the boundary layer is attached. 
This ensures the absence of any upstream influence between different regions and 
hence eliminates the otherwise costly interfacing between neighboring regions. Seg- 
mentation has been implemented in the current application to reduce the size of 
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the time dependent computational domain. Moreover, segmentation allows the use 
of the more efficient space marching procedure in any region not containing an 
embedded subsonic or axially separated flow. 

RESULTS 

The current numerical technique has been used to simulate the flow field around 
the entire configuration of the shuttle orbiter and a simple projectile with a canard. 
Shuttle Orbiter 

The supersonic viscous flow about the shuttle orbiter has been computed using 
the four computational segments or regions shown in Fig. (2). The bow shock 
segment, Region I, is treated as a typical blunt body problem with supersonic 
outflow and is solved using a time-relaxation code. Region II is centered about 
the canopy and is also solved with the time dependent code because the canopy 
shock induces a small region of subsonic flow. Region III is purely supersonic and 
could be solved by space marching, but for ease of application the more robust 
time dependent code was again utilized. Finally, Region IV containing the wing, 
the OMS pods, and vertical tail has embedded regions of subsonic flow as well as 
regions of axial separation and must be solved using the time dependent code. In 
Region IV it is particularly important that the generated mesh be body normal, 
both in the cross-sectional and axial planes in order to resolve the viscous flow over 
different components of the orbiter which are inclined at relatively steep angles to 
the fuselage. 

The hyperbolic grid generation procedure was used to generate the grids in Re- 
gions II, III and IV. whereas an algebraic procedure was used for the simple nose 
geometry. Figure (3) indicates the 25 > 123 x 29 grid generated for Region IV in 
which there are 25 points in the axial direction, 123 points in the circumferential 
direction, and 29 points in the normal direction. Projections of either r) — constant 
or onto £ = constant grid planes onto either y = constant or x = constant planes 
are used to illustrate the grid. Figure (3a), for example, shows a segment of the 
grid near the beginning of the vertical tail in the r\ = constant leeward plane of 
symmetry. Figures (3b) to (3d) show various £ = constant grid planes projected 
into the y - z plane. It is seen that the grid lines are normal to the body surface. 

The pressure contours in a cross section containing the wing, OMS pods and 
vertical tail are shown in Fig. (4). A side view of the simulated oil flow (i.e. 
limiting streamlines) at the body surface are shown in Fig. (5). It is seen that there 
is an area of axially separated flow upstream of the vertical tail and the OMS pods. 
The location of the separation line coincides with a line along which the limiting 
streamlines coalesce while the reattachment line coincides with a line at which the 
limiting streamlines diverge. 

The computations were performed on the CRAY-XMP-12 with one central pro- 
cessor and two millions words of memory. The code is vectorized and the metrics 
are not stored but are computed when needed. The computer time required per 
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grid point is 0.00012 sec per iteration. 

Canard Projectile 

The numerical procedure was also used to compute the fiowfield around a projec- 
tile with a canard at a Mach number of 4 and zero angle of attack. The geometry 
(Fig. (6)) consists of an 8 degree blunted cone followed by a cylinder. A canard 
protrudes from the cylinder at a plane 90 degrees from the pitch plane of symmetry. 
The hyperbolic grid generation scheme was used to create the three- dimensional 
mesh required for the viscous computations. Various views of the grid are depicted 
in Fig. (7). It is seen that the grid is body normal in both the streamwise and 
meridional directions. It is also seen that a clustering is employed in the outward 
direction in order to resolve the viscous effects near the body surface. In addition, 
a simple streamwise clustering is introduced to adequately resolve the flow features 
in the wake of the canard. 

In generating the grid, the cell volumes were specified such that the outer bound- 
ary is nearly conical with a half cone angle of approximately 20 degrees. This assures 
that the bow shock is located well inside the outer boundary because the shock angle 
corresponding to an 8 degree cone at Mach number 4 is about 18 degrees. 

The above mesh was then partitioned into three longitudinal segments. The 
first segment contains the spherical nose. The second segment contains the conical 
section and the front portion of the canard. The third segment contains the back 
end of the canard and the cylindrical portion. An axisymmetric blunt body code 
was used to compute the flow field in the first segment. The flow field in the second 
and third segments were obtained using the time relaxation procedure of Ref. ;3;. 
Because the body is at zero angle of attack, the computational domain was confined 
to the region between 4 >= 0 and 6 = 90 degrees. The results of the calculations in 
the second segment revealed the absence of any embedded subsonic flow bubbles 
or axially separated flow regions. This means that a PNS code could have been 
utilized in this region. 

The results of the computations in the third segment are shown in Figs. (8-12). 
The Mach contours in the q>= 90 degrees symmetry plane (canard plane) are shown 
in Fig. (8). It is seen that there is a relatively large subsonic bubble residing in 
the canard wake. The velocity vectors in the canard plane are shown in Fig. (9). 
The regions of axial flow separation and reattachment are clearly visible from the 
corresponding limiting streamlines depicted in Fig. (10). The pressure contours and 
velocity vectors in a streamwise stations corresponding to the canard full extent and 
in the cylindrical portion are shown in Figs. (11, 12) and (13, 14) respectively. The 
above computations were obtained using 88125 grid points on the CRAY-XMP- 
12. The computer time per time step is approximately 13 seconds. A total of 800 
iterations were required for convergence. 

Summary 

The present paper described a numerical mesh generation technique to be used 
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with an implicit finite-difference method for simulating viscous supersonic flow 
about low-aspect-ratio wing body configurations using a single grid strategy. The 
computational domain was segmented into multiple regions, with borders located in 
supersonic areas to avoid the otherwise costly interfacing procedure between adja- 
cent segments. The numerical procedure has been successfully applied to calculate 
the turbulent flow around the shuttle orbiter and a canard projectile at supersonic 
free-stream Mach number. 
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Fig. 3c Projection of the grid on an axis normal Fig. 3d Projection of the grid on an axis normal 
plane containing the wing. plane containing the wing. OMS pods and vertical 
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Fig. 4 Pressure contours in a cross section near Fig. 5 Side view of the oil flow (li 
the end of region IV. lines) in region IV. 
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Fig. 6 Flow regions around a projectile with a Fig. 7a Grid in the canard plane of symmetry, 
canard. 



Fig. 7b Grid in the leeward plane of symmetry. Fig. 7c Grid in a cross section containing the 

canard. 
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Fig. 8 Mach number contours in the canard sym- Fig. 9 Velocity vectors in the canard plane of 
metry plane. symmetry. 



Fig. 10 Limiting streamlines (simulated oil flow) 
in the canard region. 
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Fig. 11 Pressure contours in a cross section cor- 
responding to the full extent of the canard. 


Fig. 12 Velocity vectors in a cross section corre 
sponding to the full extent of the canard. 


Fig. 13 Pressure contours in a cross section in 
the wake of the canard. 


Fig. 14 Velocity vectors in a cross section in the 
wake of the canard. 
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